%==================================================================================================
% Matlab code for tables in online Appendix F.3 in
% "Marital matching, economies of scale and intrahousehold allocations" 
% by Laurens Cherchye, Bram De Rock, Khushboo Surana and Frederic Vermeulen. 
%==================================================================================================

addpath(genpath('C:\Users\u0100562\Downloads\YALMIP-master'));
% replace the directory with the path where YALMIP is downloaded

%% loading the data
clear all
load('psid_cleaned_onlycouples.mat')

%% Bounds on economies of scale and RICEB measures
clc
ub = 0.05;               % upper bound on degree of publicness of leisure
DivorceCosts =[];     % divorce costs
bounds = [];          % economies of scale
ricebcomp = [];       % RICEB components
a = [];               % A matrix  

while ub < 0.2

    temp_id = id;

while isempty(temp_id)==0  
    fprintf('Working on marriage market id %d with Aleisure upper bound as % d\n',temp_id(1), ub);
    %preparing the indices which are in the marriage market being evaluated
    id_ind = find(id==temp_id(1));
    display('Computing Divorce Costs');
    [temp_DivorceCostsM,temp_DivorceCostsF,temp_DivorceCostsMF] = get_divorce_costs(leisurem(id_ind),leisuref(id_ind), hworkm(id_ind),hworkf(id_ind), q(id_ind),wagem(id_ind),wagef(id_ind),nonlabor(id_ind),ub);
    display('Computing bounds');
    [temp_min_scale,temp_max_scale,temp_mina,temp_maxa,temp_min_scalem,temp_min_privatescalem,temp_min_publicscalem,temp_max_scalem,temp_max_privatescalem,temp_max_publicscalem,temp_min_scalef,temp_min_privatescalef,temp_min_publicscalef,temp_max_scalef,temp_max_privatescalef,temp_max_publicscalef] = get_scale(leisurem(id_ind),leisuref(id_ind), hworkm(id_ind), hworkf(id_ind), q(id_ind),wagem(id_ind),wagef(id_ind),nonlabor(id_ind),temp_DivorceCostsM,temp_DivorceCostsF,temp_DivorceCostsMF,ub);
   
    % merging bounds on economies of scale  and RICEB measures
    bounds_temp = zeros(length(id_ind),11);
    bounds_temp(:,1) = temp_min_scale';                        % minimum economies of scale
    bounds_temp(:,2) = temp_max_scale';                        % maximum economies of scale
    bounds_temp(:,3) = temp_max_scale' - temp_min_scale';      % difference in max and min economies of scale
    bounds_temp(:,4) = temp_min_scalem';                       % minimum male RICEB 
    bounds_temp(:,5) = temp_max_scalem';                       % maximum male RICEB
    bounds_temp(:,6) = temp_max_scalem' - temp_min_scalem';    % difference in max and min male RICEB  
    bounds_temp(:,7) = temp_min_scalef';                       % minimum female RICEB
    bounds_temp(:,8) = temp_max_scalef';                       % maximum female RICEB
    bounds_temp(:,9) = temp_max_scalef' - temp_min_scalef';    % difference in max and min female RICEB
    bounds_temp(:,10) = familyid13(id_ind);                    % id of the household 
    bounds_temp(:,11) = ub*ones(length(id_ind),1);             % upper bound on Aleisure used 
    bounds = [bounds; bounds_temp];                         
    
    % merging RICEB components
    ricebcomp_temp = zeros(length(id_ind),11);                  
    ricebcomp_temp(:,1) = temp_min_privatescalem';             % share of private consumption when male RICEB is minimized
    ricebcomp_temp(:,2) = temp_min_publicscalem';              % share of public consumption when male RICEB is minimized
    ricebcomp_temp(:,3) = temp_max_privatescalem';             % share of private consumption when male RICEB is maximized
    ricebcomp_temp(:,4) = temp_max_publicscalem';              % share of public consumption when female RICEB is maximized
    ricebcomp_temp(:,5) = temp_min_privatescalef';             % share of private consumption when female RICEB is minimized
    ricebcomp_temp(:,6) = temp_min_publicscalef';              % share of public consumption when female RICEB is minimized
    ricebcomp_temp(:,7) = temp_max_privatescalef';             % share of private consumption when female RICEB is maximized
    ricebcomp_temp(:,8) = temp_max_publicscalef';              % share of public consumption when female RICEB is maximized
    ricebcomp_temp(:,9) = familyid13(id_ind);                  % id of the household
    ricebcomp_temp(:,10) = ub*ones(length(id_ind),1);          % upper bound on Aleisure used 
    ricebcomp = [ricebcomp; ricebcomp_temp];
    
    % merging the elements of A matrix 
    a_temp = zeros(length(id_ind),12);                     
    a_temp(:,1:5) = temp_mina;                                 % components of A matrix when economies of scale is minimized
    a_temp(:,6:10) = temp_maxa;                                % components of A matrix when economies of scale is maximized
    a_temp(:,11) = familyid13(id_ind);                         % id of the household
    a_temp(:,12) = ub*ones(length(id_ind),1);                  % upper bound on Aleisure used 
    a = [a; a_temp]; 
    
    % merging divorce costs 
    DivorceCosts_minm = min(temp_DivorceCostsMF,[],2);          
    DivorceCosts_minf = min(temp_DivorceCostsMF)';             
    DivorceCosts_avg = (nansum(temp_DivorceCostsMF)' + nansum(temp_DivorceCostsMF,2))./(sum(~isnan(temp_DivorceCostsMF))' + sum(~isnan(temp_DivorceCostsMF),2));
    
    DivorceCosts_temp = zeros(length(id_ind),8);                    
    DivorceCosts_temp(:,1) = DivorceCosts_minm;                  % minimum divorce cost for male's no blocking pairs condition 
    DivorceCosts_temp(:,2) = DivorceCosts_minf;                  % minimum divorce cost for female's no blocking pairs condition  
    DivorceCosts_temp(:,3) = min(DivorceCosts_temp(:,1:2),[],2); % minimum divorce cost for couple's overall no blocking pairs condition 
    DivorceCosts_temp(:,4) = DivorceCosts_avg;                   % average divorce cost for couple's overall no blocking pairs condition
    DivorceCosts_temp(:,5) = temp_DivorceCostsM;                 % divorce cost for male's individual rationality condition 
    DivorceCosts_temp(:,6) = temp_DivorceCostsF;                 % divorce cost for female's individual rationality condition
    DivorceCosts_temp(:,7) = familyid13(id_ind);                 % id of the household
    DivorceCosts_temp(:,8) = ub*ones(length(id_ind),1);          % upper bound on Aleisure used 
    DivorceCosts = [DivorceCosts; DivorceCosts_temp];               
    
    % go to the data of next marriage market
    temp_id(find(temp_id==temp_id(1)))=[];                       % remove these households and go the households in next marriage market
end

fprintf('Computations for ub equal to %d are finished \n', ub);
ub = ub + 0.05;
end


minscale_5 = bounds(1:1322,1);
maxscale_5 = bounds(1:1322,2);
diffscale_5 = bounds(1:1322,3);
minscalem_5 = bounds(1:1322,4);
maxscalem_5 = bounds(1:1322,5);
diffscalem_5 = bounds(1:1322,6);
minscalef_5 = bounds(1:1322,7);
maxscalef_5 = bounds(1:1322,8);
diffscalef_5 = bounds(1:1322,9);

minscale_10 = bounds(1323:2644,1);
maxscale_10 = bounds(1323:2644,2);
diffscale_10 = bounds(1323:2644,3);
minscalem_10 = bounds(1323:2644,4);
maxscalem_10 = bounds(1323:2644,5);
diffscalem_10 = bounds(1323:2644,6);
minscalef_10 = bounds(1323:2644,7);
maxscalef_10 = bounds(1323:2644,8);
diffscalef_10 = bounds(1323:2644,9);

minscale_15 = bounds(2645:3966,1);
maxscale_15 = bounds(2645:3966,2);
diffscale_15 = bounds(2645:3966,3);
minscalem_15 = bounds(2645:3966,4);
maxscalem_15 = bounds(2645:3966,5);
diffscalem_15 = bounds(2645:3966,6);
minscalef_15 = bounds(2645:3966,7);
maxscalef_15 = bounds(2645:3966,8);
diffscalef_15 = bounds(2645:3966,9);

save psid_results_f3_onlycouples;




